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ABSTRACT 

We discuss the basic hydrodynamics that determines the density structure of 
the disks around hot stars. Observational evidence supports the idea that these 
disks are Keplerian (rotationally supported) gaseous disks. A popular scenario 
in the literature, which naturally leads to the formation of Keplerian disks, is the 
viscous decretion model. According to this scenario, the disks are hydrostatically 
supported in the vertical direction, while the radial structure is governed by the 
viscous transport. This suggests that the temperature is one primary factor that 
governs the disk density structure. In a previous study we demonstrated, using 
3-D NLTE Monte Carlo simulations, that viscous keplerian disks can be highly 
non-isothermal. In this paper we build upon our previous work and solve the full 
problem of the steady-state non- isothermal viscous diffusion and vertical hydro- 
static equilibrium. We find that the self-consistent solution departs significantly 
from the analytic isothermal density, with potentially large effects on the emer- 
gent spectrum. This implies that non-isothermal disk models must be used for a 
detailed modeling of Be star disks. 

Subject headings: radiative transfer — stars: emission line, Be — circumstellar 
matter — polarization 
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Introduction 



Disks are observed in a wide variety of astrophysical systems, and they are created in 
both inflows and outflows whenever the rotation speed greatly exceeds the radial flow speed. 
Examples of disk systems formed during infall include: young stellar objects, mass transfer 
binaries, and active galactic nuclei. Similarly, disks also form in the outflowing winds of 
luminous stars, such as Be stars, B[e] stars, and possibly LBV and AGB stars. From the 
kinematic point of view, the disks of hot stars are probably relatively simple, so they form 
ideal laboratories for studying the ph ysics of outflowing circum stellar disks. For a general 
review of Be stars and their disks, see iPorter fc Riviniud (120031 ). 



The current paradigm for Be star di sks is that of a geometrically th in viscous Keplerian 



disk in vertical hydrostatic equilibrium (ILee et alj 



been corroborated by several recent studies (e.g. ICarciofi et al 



1991 



Qkaz aki 200 l|). This paradigm ha s 



20061 : iMeilland et alJl2007h . 



The time-independent solution for the viscous decretion disk str ucture is relatively straight- 



forward to ob tain if one assumes the disk to be isothermal (e.g. iBjorkman & Carciofil 12005 



Okazakill200ll ). In this case, one finds that the disk is in vertical hydrostatic equilibrium and 
because the rotation velocity is much larger than the sound speed, the disk is geometrically 
very thin, with opening angles of a few degrees only. In addition, the solution of the viscous 
outflow results in a simple power-law for the radial dependence of the gas density. 



In a previous paper (ICarciofi &: Bjorkmanll2006l . hereafter Paper I), we used the three- 
dimensional (3-D), non-local thermodynamic equilibrium (NLTE) code HDUST to study 
the temperature structure of Keplerian disks. We found that those disks can be highly non- 
isothermal, mainly in their denser inner parts. Later, our results were largely confirmed by 
Sigut fc Jonesl (120071 ). who studied the thermal properties of the Be star 7 Cas. 



In Paper I we adopted the isothermal analytic density structure for the disk. This is 
clearly inconsistent with the fact that the disks are non- isothermal, and it is the purpose of 
this paper to remove this inconsistency. We present the first self-consistent solutions for the 
structure of viscous decretion disks, taking into full account all the relevant physical processes 
involved. This is a somewhat complicated problem, because the disk temperature controls 
the geometry (via the hydrostatic equilibrium and viscous diffus ion equations), which in tur n 
determine the heating and, hence, the temperature of the disk (IKenyon Hartmannlll987l ). 
Solving this problem only became feasible after the advent of large computer clusters. 

We also study the emergent spectra of representative model and compare them with 
those of isothermal models to understand to what extent the non-isothermal solution is 
required for fitting actual observations. 
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2. Disk Theory 



The v iscous decretion model is essentially the same as that employed for protostellar 
disks (e.g.. lBjorkmanlll997l ). the primary difference being that Be disks are outflowing, while 



pre-main-sequence disks are inflowing. 

It is supposed that some (as yet unidentified) mechanism injects material with Kep- 
lerian velocities into the base of the disk. Provided this material injection is continual, 
eddy/turbulent viscosity then transports angular momentum from the inner boundary of 
the disk outwards. 

The derivation of the density structure for isothermal disks is well -known. Below, we 
solve t his problem for the non-isothermal case. We follow the derivation of lBjorkman fc Carciofi 
( 120051 ). assuming that the temperature is axisymmetric, i.e., T = T(w,z). 



2.1. Hydrostatic Structure 

Let us, initially, ignore viscosity. In cylindrical coordinates (w,4>,z), the steady-state 
fluid equations are 

Id. . 1 d , , d . . . . 

{wpv^)-\ -kt{pvj>) + -^-{pVz) = 0, (1) 



w dw w d<p dz 

dv^ dv m dv m v\ ldP 

Vw-~— + aT+ v *-a = a - + f™ ' 2 

ow w ocp oz w p aw 

dvtf, v^dv^ dv^ v^v^ 1 OP 

Vvjjt- + — -E~ + V z^r + — - = -ai + U > 3 

aw w ocp oz w pw ocp 

dv z v (f) dv z dv z ldP 

V *>^- + —^T + V *^- = — -rr + fz- (4) 
aw w ocp oz p oz 

where p is the gas density, P is the pressure and v w , and v z are the velocity components. 
In the absence of any other force, the external force components are given by the gravity of 
the central star 

GMw 

J* (^2 + ^2)3/2 ' 

GMz 

Jz {w^ + z^ ' K ' 

Assuming circular orbits {v w = 0, and v z = 0), the only non-trivial fluid equations are 

p ow w 
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IdP 

p dz 



fz 



the w- and z-momentum equations, respectively. To specify the pressure, we introduce the 
equation of state P = a 2 p , where a = a(w, z) is the sound speed that is related to the gas 
kinetic temperature by a = (kT) 1 ^ 2 (pm-a)" 1 ^ 2 . In this last expression, k is the Boltzmann 
constant, p is the gas molecular weight and mn is the mass of the hydrogen atom. . 



In the thin-disk limit (z <C w), we obtain 



d\n(a 2 p) 



Krit {R/ruf 12 



V^Rz 
a 2 w 3 



(9) 
(10) 



dz 

1 /2 

where the critical velocity, V CT it = (GM/R) , is the Keplerian orbital speed at the stellar 
surface and R is the stellar radius. 

The above equations mean that the disk rotates at the Keplerian orbital speed and is 
hydrostatically supported in the vertical direction. Defining the disk scaleheight 



H {W) = W ' 

t'crit 



(11) 



we obtain, by integrating the hydrostatic equilibrium equation ( TlOT) along the ^-direction 
and writing the result in terms of the temperature 



In 



T(zu, z)p(tu, z) 
T (w)po(w) 



T(w,z') 



where the subscript indicates values for z — 0. Solving for p(zu, z) we obtain 



p(w,z) = pa(w) 



T (zu) 
T(w,z] 



exp 



2H 2 (w) 



T(zo,z') 



(12) 



(13) 



As we shall see below, it is useful to express p(w, z) in terms of the disk surface density, 
H{w), written as 

/oo 
p{vj,z)dz . (14) 
-oo 



Substituting eq. f[T31 into eq. (IT4T) and solving for po(zu) we obtain 



p(w,z) 



E(ot) 
T(zu, z) 



exp 



UW) dz' 2 



2H 2 (w)J T(m,z> 



T(w,z] 



exp 



T (U7) 2 



2H 2 (w)J T{w,z>) 



dz 



(15) 
(16) 
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2.2. Viscous Outflow 

Clearly, in eq. ffT6l) the disk density scale £(W) is completely undetermined, because 
for a pure Keplerian disk we can choose to put an arbitrary amount of material at a given 
radius, w. To determine the density scale, we must include a mechanism that transports 
material from the star outwards. The most popular way to generate decretion is to consider 
viscous effects. 

If we add viscosity to the fluid equations ([I]) to (j4]), we can still assume that the disk is 
axisymmetric and that the vertical structure is hydrostatic (v z = 0). However, the presence 
of an outflow implies that v w ^ 0. The w- and z- momentum equations are the same as 
before, so and p are the same as in the pure Keplerian case [eqs. ([9]) and (JTB]) ]. 

The continuity equation (TTJ becomes 

= , (17) 

ozu 

so the decretion rate M = 2itvjY>v. D j is a constant (independent of w). The viscous outflow 
speed is given by 

M , . 

v„ = t-^; ■ (18) 

The remaining fluid equation is the 0-momentum equation ([3]). This equation now is more 
complicated because viscosity exerts a torque, which is described by the viscous shear stress 
tensor, 7Tjj. Including this shear stress, the 0-momentum equation becomes 

dvj, v^v^ 1 <9 2 
V -d^ + — = ^d^ {W7T ^ ' (19) 

where 

= ,P«%^ . (20) 



For th e kinematic viscosity, z/, we use the prescription introduced by IShakura fc Sunyaev 



(119731 ). according to which v = aaH, where < a < 1. The viscous shear stress tensor then 
becomes 

^ = ~-aa 2 p . (21) 
Multiplying eq. ( fl9l) by pw 2 and integrating over and z, we find 

■ d d f°° 

—M— — \wv<f) = — — / G7vr ro 2ttw dz . (22) 



dzu dm 

Substituting eq. fTST]) into eq. fl22|) and integrating it over w we obtain 



' p{m, z)dz + constant , (23) 
-oo T (w) 
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where the integration constant is determined by the outer boundary condition of the disk. 
Finally, defining the weighted mass average temperature 

/rrf U _ IZ T{zu,z)p{w,z)dz J™ oo T(vu,z)p(vu,z)dz 

{T{w)) = — — = —— , (24) 

substituting eq. (124")) into eq. (|23|) . and solving for the surface density we find 

E(U7) = tok«(T{m))w*l* [{RdM 1] ' (25) 

where po is the mid-plane mean molecular weight, and the disk outer radius, Rd, is deter- 
mined by the integration constant in eq. ( I23i) . Physically this arises from the outer boundary 
condition that sets the sonic location (critical point) in the disk outflow. If the disk is trun- 
cated by a physical mechanism, such as tidal truncation by a binary, this location will be 
near the truncation radius of the disk. 

We have now completed the hydrodynamic description of a non-isothermal viscous de- 
cretion disk. To fully determine the problem, one must specify the decretion rate, M, the 
disk radius, Rd, and the disk temperature structure, T{w,z). 



2.3. Properties of the Isothermal Solution 



The analytical solution for the isothermal case is fully discussed elsewhere (e.g. iBjorkman fc Carcion 



20051 ; IOkazakill200ll ). but since we will focus below on the comparison between the isothermal 



and non-isothermal cases, it is useful to present a brief outline of its properties. 

Assuming a constant temperature we find, from eq. (125]) . that the surface density de- 
pends on radius as 

S(w) oc zu- 3 / 2 [(R A /w) l/2 - 1] . (26) 
For large disks (R -C -Rd), the surface density becomes a simple power-law, T>(zu) oc uj~ 2 . 

The isothermal disk flares quite dramatically with radius, H oc zu 1 ' 5 . Since p oc T>/H, 
we obtain that the isothermal disk density profile is quite steep, p oc ro~ 3 - 5 . 

Finally, we find that, for large disks, the radial velocity is a linear function of the radial 
distance, v m oc tu. 
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3. Self-consistent Hydrostatic Solutions 



In this section we determine numerically the solution of non-isothermal viscous decretion 
disks, and study how it compares with the simple analytical solutions of the isothermal case. 
For this, we modified the computer code HDUST described in Paper I. 

HDUST is a fully 3-D NLTE Monte Carlo radiation transfer code that solves the radia- 
tive equilibrium temperature for arbitrary density and velocity configurations. In brief, the 
NLTE Monte Carlo simulation performs a full spectral synthesis by emitting stellar photons 
with random f requencies sampled using a Kurucz model atmosphere for the central star 
( iKurucz Ill994l ). Each photon is followed as it travels through the envelope (where it may 
be scattered, or absorbed and reemitted, many times) until it escapes. 

During the simulation, whenever a photon scatters, it changes direction, Doppler shifts, 
and becomes partially polarized. Similarly, whenever a photon is absorbed, it is not de- 
stroyed; it is reemitted locally with a new frequency and direction determined by the local 
emissivity, j u , of the gas. Note that we include both continuum processes and spectral lines 
in the opacity and emissivity of the gas. Since photons are never destroyed (absorption is al- 
ways followed by reemission of an equal energy photon packet ), our procedure automati cally 
enforces radiative equilibrium and conserves flux exactly (see iBjorkman fc Woodll200ll ). 



The interaction (absorption) of the photons with the gas provides a direct sampling of 
all the radiative rates, as well as the heating and cooling of the free electrons. Consequently, 
we adopt an iterative scheme in which at the end of each iteration we solve the rate equations 
to update the level populations and electron temperature of the gas. The iteration proceeds 
until convergence of all state variables is reached. 

In paper I we discussed in great detail the solution for the temperature structure of a 
viscous decretion disk with a fixed density structure. We found that the disk temperature 
is a hybrid between that of optically thick disks of Young Stellar Objects and optically thin 
winds from Hot Stars. For the optically thin upper layers, as well as for the outer disk 
(w/R> 10) , the temperature is nearly isothermal, with an average temperature of about 
60% of T eff . 

At the mid-plane the temperature has a complex structure. Close to the photosphere the 
temperature is a yery steep function of tu, with a profile well described by a flat blackbody 
reprocessing disk (j Adams et al.lll987l ). 

1/4 



T d (w) 



71 



1/4 



R 

sm *| — 

V3 



R 



R 2 



W V 



(27) 



where T* is the temperature of the radiation that illuminates the disk. When the disk 
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becomes optically thin vertically, the temperature departs from this curve and rises back to 
the optically thin radiative equilibrium temperature. 

For the present work we modified HDUST to solve simultaneously for both the disk 
density structure and the gas state variables. At the end of each iteration, we first determine 
the new temperature and, then, solve the hydrostatic equilibrium equation, eq. (fTUj) . for the 
vertical structure of the disk, as well as the steady outflow, eq. (l25j) . for surface density. The 
new density thus obtained is used for the next iteration. The procedure is repeated until the 
density, temperature, and level populations converge (this typically requires between 5 to 30 
iterations, depending on the disk density: the larger the density, the larger the number of 
iterations). 

The primary parameter determining the disk structure is the density scale, which is 
set by the ratio of the stellar mass loss rate and the parameter a [Eq. (1251) ]. The other 
parameters, such as stellar effective temperature, stellar radius, etc., albeit important, do 
not change significantly the properties of the temperature solution (Paper I). For this reason, 
we choose a fixed set of parameters to define the star and the disk size and leave only the 
density scale as a free-parameter in this study. We define 3 density regimes: low mass loss 
rate, high mass loss rate and an equal mass case, whose meaning will be explained below. 
All fixed model parameters are listed in table [H and the parameters associated with each 
density regime are listed in table [2j 

For the parameter a we arbitrarily choose a value of 0.1 so that the viscous diffusion 



timescale roughly agrees with the observed disk dissipation and creation timescales (jOkazaki 



20011 ). However, the actual value of a is, within a reasonable range, largely irrelevant for the 
present problem. This becomes evident when one analyses the role of a in the disk dynamics. 
From Eqs. (ITSl) and (12"5|) we find that v m oc a; hence, a controls the viscous outflow speed. 
However, as long as the outflow is highly subsonic, which is true for values of a in the range 
0.1 — 1, for instance, the disk dynamics are entirely dominated by rotation. It follows that 
if we vary a within the above range while adjusting M to maintain a constant value of M/a 
(so that the disk density does not change), all the properties of the solution will be essentially 
unaltered. 

Our solutions for the disk structure are shown in Figures [1] to HI In Figure [1] we compare 
the density structure of the consistent numerical solution with the structure of a correspond- 
ing isothermal disk. For the later model, we adopted a value for the gas kinetic temperature 
of 60% of T e g, or 12 000 K (see Paper I). Two things are evident from the figure: the 
consistent solution has smaller scaleheights and larger densities in the mid-plane than the 
isothermal density. Those features are easier to see in Figure |SJ in which we plot the ratio 
between the self-consistent and isothermal solutions, p/pi SO - For the high mass loss rate 
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case (bottom plots) the differences in the inner part of the disk are quite significant: the 
equatorial densities are enhanced by a factor of up to ~ 3, whereas the upper layers of the 
disk have much lower densities. For low mass loss rates there is a slight enhancement in the 
equatorial regions, but the differences are much smaller than the previous case. 

These results can be understood in terms of the temperature structure. In the last 
section we demonstrated that what controls the steady flow solution is (T), the weighted 
mass average of the temperature, eq. ([21]) • This quantity is shown in Figure [3] for both 
the low and high mass loss rate models. For the later, (T) is large at the base of the star 
(« 16 000 K), drops very quickly to values close to 6 000 K around w/R = 5 and reaches 
a near-constant plateau of ~ 12 000 K around w/R = 10. The low values of (T) cause 
the material to collapse towards the equator, thus increasing the mid-plane densities and 
decreasing the densities of the upper layers. Note that (T) is very similar to the mid-plane 
temperature (left plots of Figure [3]), because the mid-plane region dominates the integral in 
eq. (J2U). 

For low densities there are also important variations in (T), but they are confined to 
much smaller radii and are of smaller amplitude. For this reason, the consistent solution 
does not differ significantly from the corresponding isothermal model for low mass loss rates, 
as shown in Figure [2j 

In Figure H] we exhibit other aspects of the self-consistent solution for the disk density. 
We define the disk opening angle, 9, as the angle subtented by one scaleheight at a given 
radius, 

w 

Let us first analyze the opening angle of the high mass loss rate model (Fig. HJ top right 
panel). Initially, at w/R = 1, the scaleheight is larger than the isothermal value because 
of the large temperature at the base of the disk. As the temperature quickly falls, the disk 
flares negatively, i.e., its geometrical thickness becomes smaller with growing radius. This 
is a subtle effect but quite interesting in view of the striking contrast with the isothermal 
solution. The disk begins to flare again around w/R = 2 as a result of the decrease of 
gravitational potential, but it does not flare nearly as much as an isothermal disk. When 
(T) reaches its minimum and starts rising again, the disk flares dramatically, reaching the 
isothermal opening angle when (T) becomes approximately constant. 

These remarkable properties of the disk opening angle are also present in the low density 
model, but are much less pronounced. En passant, we note that the opening angle for the 
inner disk is very small for all models, in agreement with the general belief that the disks of 
Be stars are geometrically very thin. 



w i = sin 



(28) 
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The remaining plots of Figure H] illustrate different aspects of the radial density structure. 
The plots second from bottom, for instance, show in a more quantitative way the density 
enhancement of the mid-plane, which can be larger than 2.5. In those plots, a feature that 
is hard to see in Figure [2] is quite evident: very close to the star there is a density decrement 
in the mid-plane, which arises as a result of the large temperatures there. As before, the 
effects are less conspicuous for low mass loss rates. 

Another striking feature of the solution can be seen in the bottom panels of Figure HJ 
where we plot dln(p )/dln(w) as a function of w. This quantity gives the local index of the 
density profile, were it to be fitted by a power-law. The radial dependence of the density 
departs very significantly from a simple power-law. Close to the star {w/R < 3) the density 
profile is much flatter than the isothermal profile, but after the (T) turnover point, which for 
the high mass loss rate model is around w/R = 3, the profile becomes significantly steeper. 
The density profile only approaches the isothermal profile after (T) reaches the near constant 
plateau at w/R ~ 10. 

We end this section discussing a subtle property of the solution for the disk structure. 
For a given M, the total disk mass of the consistent solution is always larger than the disk 
mass of an isothermal disk with T e = 60% T e g. For instance, for the high mass loss rate case 
(M = 5 x 10~ n M & yr _1 ) the disk mass is 3.6 x 10 -9 M Q whereas the isothermal disk mass 
is 2.8 x 10 -9 M Q . The total disk mass scales with the global weighted mass average of the 
temperature, which is lower than the assumed temperature for the isothermal model. For 
this reason, we defined the equal mass cases in Table [2] so that the total disk mass of the 
consistent numerical model equals that of the isothermal model. 

4. Observational Signatures 

In the previous section we have shown that the self-consistent solution for the disk 
structure can depart significantly from the analytical isothermal solution. Here we study 
what are the effects of the density changes on the emergent spectrum. For each density 
regime we compare models with the same total disk mass. 

To clearly separate between different physical effects, i.e., what is the effect of the density 
and what is the effect of the temperature, we define three different models. The first model, 
which we call consistent model, corresponds to the full non-isothermal solution for both the 
hydrostatic equilibrium and radial outflow. The second model, which we call mixed model, 
has a fixed density structure given by the analytical isothermal solution, but the gas state 
variables (temperature and level populations) are determined by HDUST. Those are the 
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same type of models presented in Paper I. Finally, we define isothermal models, for which 
both the density and temperature are fixed, and only the level populations are calculated. 

The comparison between the mixed and isothermal models allows us to pinpoint the 
effects of the temperature structure on the emergent spectrum; on the other hand, the 
comparison between the consistent and mixed models enables us to determine the effects of 
the density changes. 

Based on the previous description of the properties of the solution, we expect the largest 
differences between the three models to occur for higher mass loss rates (higher densities). 
Indeed, Figures [5] and [6] confirm this expectation, since the results for the low mass loss 
rate models are nearly identical, whereas the results for higher densities can be substantially 
different, mainly for the emission line strengths and the continuum polarization. 

It is important to emphasize that, for the low mass loss rate regime, the results for 
the isothermal model only match those of the mixed and consistent models because the 
proper choice of the temperature was made: 12 000 K or 60% of T cS . If another temperature 
choice had been made, the differences would be significant. We can conclude, therefore, that 
models with low mass loss rates and with the correct choice of the temperature do a very 
good job in determining the emergent spectrum of viscous decretion disks. This means that 
employing the more complex mixed or consistent solutions is not necessary (of course, the 
NLTE treatment of the level populations must still be done). 

The situation for the high mass loss rate case is the opposite, since relevant differences 
exist in the emergent spectrum of the three models. The key for understanding those dif- 
ferences resides, in part, in the Hydrogen level populations. As shown in Figure [7J at the 
non-isothermal part of the disk (m/R < 10) the level populations can be different by factors 
of up to 2 — 3. This indicates that both the temperature structure and the density structure 
affects the populations in a complex way. 

Let us first compare the emergent spectra of the mixed and isothermal models, which, 
as stated above, can give us clues on how the temperature structure affects the results. We 
recall that those models have exactly the same density distribution. There are significant 
differences in the continuum polarization (factor of about 2) and in the Ha emission (about 
50%), and those differences come directly from the dissimilar excitation structure in the 
non-isothermal part of the disk. This proves the importance of using the correct solution 
for the temperature and a full NLTE treatment of the populations. In addition, there are 
small differences in the IR spectrum: the slope is slightly flatter for the isothermal model, 
which results in larger fluxes for long wavelengths. This is a direct result of the mid-plane 
temperature structure, a subtle effect which was fully discussed in Paper I. 
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The differences between the mixed and consistent models are also significant; the con- 
sistent model produces smaller polarization levels and Ha emission than the mixed model. 
Those models differ mainly in the density structure, but there are also slight differences in 
temperature and level populations as a result of the density changes. 

Another potentially useful result is illustrated in Figure [BJ where we show the Bracket 
IR spectrum for each model. Fitting the continuum spectrum with a power-law gives us 
the following spectral indices: —3.003 ± 0.001 for the consistent model; —3.167 ± 0.001 for 
the mixed model; —3.144 ± 0.001 for the isothermal model. The NIR spectral index probes 
the density structure of the inner disk. As expected, the spectral indices of the mixed and 
isothermal models are nearly identical, since they share the same density structure. The NIR 
spectrum for the consistent model, however, has a lower spectral index because of the changes 
of the radial density structure of the inner disk. This suggests the interesting prospect of 
probing the temperature structure of the inner disk by way of detailed observations of the 
NIR continuum spectrum. 



5. Conclusions 

We have investigated how the non-isothermal temperatures alter the density structure 
of Keplerian viscous disks with steady-state outflow. For low mass loss rates (low densities) 
there is relatively little affect on the disk observables because the only changes occur in the 
innermost region of the disk where it is optically thick. However for high mass loss rates, the 
disk is optically thick out to a much larger radius, so the effects are much more pronounced. 

The primary effect is that the temperature decrease causes the disk to collapse, becoming 
much thinner in the inner regions. This collapse redistributes the disk material toward the 
equator, increasing the mid-plane density by about a factor of 3 (relative to the isothermal 
density model). As a result, the disk opening angle is relatively constant (about 2-3 deg) 
out to the location of the temperature minimum. Beyond this location the temperature rise 
causes the disk to flare rapidly, and the disk opening angle approaches that of the isothermal 
model, which has opening angles of about 10° at 20 stellar radii. 

The combination of the radial temperature structure, disk scale height (opening angle), 
and viscous transport produces a complex radial dependence for the disk density. The 
equivalent power law radial density exponent varies between n = —2 in the inner disk to 
n = —5 near the temperature minimum, eventually rising back to the isothermal value 
n = —3.5 in the outer disk. 

The net result of these changes in the disk density (which alters the hydrogen level 
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populations) is that the Ha emission strength varies by about 50%, while the instrinsic 
polarization level changes by a factor of 2. The effect on the spectral energy distribuion 
(SED) is somewhat smaller. The IR continuum level changes by about 10% with slight 
changes in the spectral index (slope). 

We conclude that detailed modeling of the observations of Be star disks, such as the 
interpretation of interferometric, spectroscopic, and polarimetric observations, all require the 
use of the non-isothermal density structure. 

This work was supported by FAPESP grant 04/07707-3 (ACC) and NSF grant AST- 
0307686 to the University of Toledo (JEB). 
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Fig. 1. — Left: Solution for the disk structure of the high mass loss rate model of Table [2j 
Right: density structure of a corresponding isothermal disk. 



-16- 




2 4 6 8 051015 20 25 30 

r/R r/R 



Fig. 2. — Comparison of the hydrostatic solution with the isothermal solution. Shown is the 
ratio between the two solutions, p/ p lso . The top panels correspond to the low mass loss rate 
model and the bottom panels to the high mass loss rate models. The dashed lines correspond 
to the curves \z\ — H isot (zu). 
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Fig. 3. — Radial temperature structure. Left panels: Temperature along the mid-plane. 
Right panels: Weighted-average of T in the vertical direction, eq. (|24|) . The top panels 
correspond to the low mass loss rate model and the bottom panels to the high mass loss rate 
model. 
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Fig. 4. — Disk density structure of low (left) and high mass loss rate models (right). Top: 
Opening angle. The thick line corresponds to the self-consistent solution for the disk opening 
angle [eq. ( l28l) ]. The dashed lines correspond to opening angles for scale-heights parame- 
terized as H{w) oc zu 13 ; from bottom to top, (3 = 1.0 to 1.5, in steps of 0.1. Second from 
top: Density profile. We compare the consistent density calculated by HDUST (thick line) 
with the isothermal density (thin line). Second from bottom: Ratio between the consistent 
density and the isothermal density. Bottom: Local index of the density profile. 
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Fig. 5. — Emergent spectrum for the low mass loss rate model of Table [2] (black lines). The 
solid lines correspond to a viewing angle of 90° (disk seen edge-on) and the dotted lines to 
70°. Top panels: SED. Bottom left: polarization. Bottom right: Ha profile (a flat stellar 
profile was assumed). For comparison, we show in blue the results for the mixed model and 
in green the results for the isothermal model. For the latter models we adopted an electron 
temperature of 60% of T ef j. The red curves correspond to the stellar SED, without the disk. 
We indicate in each panel the average error of the Monte Carlo simulation, that ranges from 
0.2 — 2%. 
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Fig. 6. — Same as Fig. [5] for the equal mass model of Table [21 
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Fig. 7.— Hydrogen level populations for the high mass loss rate case. Shown are the 
level populations for n = 2 (black), and n = 3 (gray) as a function of radius. Solid lines: 
consistent model. Dotted lines: mixed model. Dashed lines: isothermal model. Each panel 
shows the populations at a given latitude (angle above mid-plane), as indicated. 



-22 - 




k 0.012 

o.oioh 



0.008 - 



0.006 - 



1.6 



1.8 2.0 
A(jLtm) 



2.2 



Fig. 8. — Bracket continuum for the high mass loss rate models of Table [2] (black line). For 
comparison, we show the mixed (dark gray) and isothermal (light gray) models. A power-law 
fit of the continuum flux is also shown as dashed lines. The average error of the Monte Carlo 
simulation is 0.3%. 
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Table 1. Fixed Model Parameters 



R 




Icrit 


a Rd 


(Re) 


(K) 


(kms -1 ) 


(R) 


7 


20 000 


400 


0.1 100 



Table 2. Model Densities 



M p (zu = R) p (zu = R) a 

(Mq yr _1 ) (g cm -3 ) (g cm -3 ) 



Low mass loss rate 


5.0 


X 


10" 


-12 


4.3 


X 


io- 


-12 


5.5 


X 


io- 


-12 


Equal Mass 


4.9 


X 


10" 


-12 


4.2 


X 


10" 


-12 


5.4 


X 


10" 


-12 


High mass loss rate 


5.0 


X 


10" 


-11 


3.0 


X 


io- 


-11 


5.5 


X 


io- 


-11 


Equal Mass 


3.9 


X 


10" 


-11 


2.5 


X 


io- 


-11 


4.3 


X 


io- 


-11 



a Value for an isothermal model 



